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1 Abstract 

Air traffic control automation synthesizes aircraft 
trajectories for the generation of advisories. 
Trajectory computation employs models of aircraft 
performances and weather conditions. In contrast, 
actual trajectories are flown in real aircraft under 
actual conditions. Since synthetic trajectories are 
used in landing scheduling and conflict probing, it 
is very important to understand the differences 
between computed trajectories and actual 
trajectories. This paper examines the effects of 
aircraft modeling errors on the accuracy of 
trajectory predictions in air traffic control 
automation. Three-dimensional point-mass aircraft 
equations of motion are assumed to be able to 
generate actual aircraft flight paths. Modeling 
errors are described as uncertain parameters or 
uncertain input functions. Pilot or autopilot 
feedback actions are expressed as equality 
constraints to satisfy control objectives. A typical 
trajectory is defined by a series of flight segments 
with different control objectives for each flight 
segment and conditions that define segment 
transitions. A constrained linearization approach is 
used to analyze trajectory differences caused by 
various modeling errors by developing a linear time 
varying system that describes the trajectory errors, 
with expressions to transfer the trajectory errors 
across moving segment transitions. A numerical 
example is presented for a complete commercial 
aircraft descent trajectory consisting of several 
flight segments. 
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2 Nomenclature and Acronyms 

A t By C Terminal sensitivity matrices 

a Speed of sound 

Ci , Cd Coefficients of lift and drag 

D Aerodynamic drag 

d Disturbances vector 

F Equations of motion 

G Constraint equations (control equations) 

hiyh p Altitude (geometric, pressure) 

i e Engine incidence angle 

L Aerodynamic lift 

M Mach number 

m Total aircraft mass 

N Output equations 

p Static air pressure 

p Parametric error vector 

S Reference wing area 

T Engine thrust 

t Time 

u Control vector 

Vt True airspeed vector 

Vi Inertial Velocity vector 

Vcas Calibrated airspeed 
V g Ground speed 

W x Easterly component of wind 

W y Northerly component of wind 

Wh Vertical component of wind 

W\ Horizontal component of wind 

W v Wind component in velocity direction 

W* Wind component in lateral velocity direction 
W y Wind component in vertical velocity direction 

£, y State and output vectors 

East and North positions 


1 



a Angle of attack 

a t Angle between thrust and relative wind 

e Parametric variations in model 

it Engine control variable (E.G. EPR, PLA, Nl) 
7 <, Flight path angle (relative to air mass) 

7 ,• Flight path angle (relative to ground) 

\p Segment termination condition 

0, V Segment transition matrices 

V a Velocity Heading (relative to air mass) 

Vi Velocity Heading (relative to ground) 

V w Wind direction 

<f> Bank Angle 

CTAS Center-TRACON Automation System 

FMS Flight Management System 

TRACON Terminal Radar Approach Control 


3 Introduction 


Increase in nation-wide air travel has put severe 
burdens on the current air traffic control system 
and controller’s workload. It has become 
increasingly important to develop computerized 
automation to assist air traffic controllers. 

Over the last decade , researchers at NASA Ames 
Research Center have designed and tested an air 
traffic control automation system, called 
Center-TRACON Automation System (CTAS). 
CTAS generates and displays descent advisories for 
air traffic controllers. Following these advisories, 
aircraft can descend efficiently, land in a scheduled 
order, and avoid potential conflict. 

In essence, an air traffic control automation 
system, such as CTAS, consists of a scheduling tool 
and a trajectory synthesizing tool. Trajectory 
synthesis computes flight trajectories for many 
aircraft. These trajectories must be consistent with 
the air traffic control regulations and pilot 
procedures. They are efficient, conflict-free, and 
meet scheduled times of arrival. Based on these 
trajectories, the scheduler generates an efficient 
landing order and landing times. Finally, the ATC 
automation system extracts advisories from these 
trajectories for each aircraft in terms of heading, 
speed, and altitude. Therefore, trajectory synthesis 
is a core element to air traffic control automation. 

Trajectory synthesis uses models of aircraft 
performance and weather conditions. These models 
necessarily contain errors from various sources. As 
a result, computed trajectories are different from 


aircraft trajectories flown in actual conditions. 
Because computed trajectories are used in conflict 
probing and scheduling, it is crucial to estimate the 
range of differences between computed trajectories 
and actual trajectories caused by modeling errors. 

This paper presents an efficient approach for 
analyzing these differences. In this paper, pilot 
feedback controls are expressed as equality 
constraints that follow ATC and pilot procedures. 
Aircraft trajectories consist of a sequence of flight 
segments defined by changing control objectives 
and segment termination conditions. Various 
sources of errors are examined and described by 
either uncertain parameters or uncertain functions. 
A constrained linearization approach is developed 
for a complete descent profile similar to profiles 
used in the CTAS system. As a result, trajectory 
differences are governed by a set of linear, 
time-varying, ordinary differential equations with 
expressions that account for the change in segment 
transition times. The proposed method offers 
physical insight into the effects of various modeling 
errors and is computationally efficient. 

In the rest of the paper, three-dimensional 
point-mass aircraft equations are first presented. 
Different components of the trajectory differences 
are analyzed, followed by discussions of 
mathematical representations by trajectory states, 
controls, and modeling errors. Then, theories of 
constrained linearization are developed. The use of 
constrained linearization is demonstrated with an 
example of a Boeing 757 descent. Time histories of 
trajectory error are demonstrated for a drag 
modeling error. 


4 Equations of Motion 


The point mass aircraft equations of motion with a 
dynamic wind environment are listed below. The 
equations of motion are written in the air-mass 
relative frame of reference, with the kinematic 
equations written in the inertial frame of reference. 


Vt 

*a 

7a 


T cos at — D 


m 


- g sin 7 a - Wy 


(T sin a t + L)sin^ W* 


mV t cos 7 a 
(T sin at + L)coe^ 
mV t 


Vt COS 7 0 
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_£COB7o —W 

V, + V, Wy 

ii = Vt cos y a sin — Wj sin 
J/i = Vt cos 7a cos a — Wt COS ^ 

A.* = Vt sin 7 a -I- W/> 

m = —my 

where, 

a t = a + t e 

W = y/wF+wj 

Other quantities of interest can be represented as 
functions of the state variables. 


M = V t /a(hi) 



Vi = (*?+*?+*?)* 


The aerodynamic and propulsive forces are 
functions of the states and the control variables 
Cl, n. 

L = ±p(A)V?SCx 

D = \p{.h)V?SC D {C L ,M) 

T = T(x,h p ,M,&ISA) 

rhf = m/(^r, A p , Af, A/5-4) 

The equations of motion form a 7 state system (Vt, 
7 a> y»» A,, m) with 3 control variables (Cl, 

it, ^). The other quantities of interest are 
represented as functions of the state variables and 
are considered to be outputs of the system of 
dynamic system. 
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All error sources that cause predicted aircraft 
trajectories to deviate from actual flown 
trajectories can be organized into the following six 
categories. 

1. Theoretical approximations 


Traiectorv Deviations 


2. Numerical approximations 

3. Modeling errors 

4. Measurement errors on the ground 

5. Measurement errors on the airplane 

6. Tracking properties of the pilot, 

autopilot, or Flight Management 

System. 

Many of the error sources can be incorporated into 
the equations of motion as additional inputs to the 
system of equations. The inputs that represent 
these errors sources are categorized by whether the 
error is constant over the course of the flight, or 
varies. The constant error sources will be collected 
in a vector called p, and the varying error sources 
in d(t). 

Figure 1 illustrates the definition of four different 
trajectories, yrs represents a trajectory computed 
by numerical integration of the TS equations, yrs 
represents the theoretical trajectory described by 
the TS equations of motion. y 0 represents a 
reference aircraft trajectory found by integrating 
the point mass equations of motion with no other 
errors introduced. Finally, y represents the actual 
aircraft trajectory described by the point mass 
equations of motion with all error sources 
introduced. 

The deviations between these 4 trajectories will be 
studied independently, with this paper focusing on 
the deviation caused by external error sources. 


1. Errors due to external error sources 

2. Errors due to different sets of equations 

3. Numerical integration errors 



Figure 1: Components of Trajectory Deviations 

For a more detailed description of error sources and 
components of trajectory deviations, see Reference 
[!]• 
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6 Trajectory Representation 

In both actual flight and in trajectory predictions 
using FMS and TS algorithms, a complete 
trajectory is divided into a series of segments. The 
segments are defined by 2 control objectives for the 
longitudinal motion, and 1 control objective for the 
lateral motion. The 2 longitudinal control 
objectives are chosen from 3 groups of objectives: 
speed (Mach number, CAS), vertical path (altitude, 
altitude rate, flight path angle), or throttle setting. 

The flight segment transitions are determined by 
segment termination conditions, or “capture 
conditions” . Typical capture conditions are 
altitude, CAS, Mach, and path distance. 

The system of equations is represented 
symbolically by collecting the states, controls, 
inputs, outputs, and error sources into vectors. 

X = [Vt,*a,7;*i,Vi, hi,m,t)' 

a - [*■ ,CL,<t>)‘ 

V = [Vcas, Vcas, M, M, ♦ t> 7,-,7<, V t , Vi, 

hp , hp , hp , h{ , V t , Vcas > r F\ 
d = [AWg,AW y ,AW h ]' 

AV e , AM e , Ah e , Axtod ] 1 

The perturbed equations of motion can be 
represented by the system of equations: 


X = 

F(x,^,d,P) 

(1) 

X(t ic ) = 

Xic + AXi e 

(2) 

0 = 

G(x,u,d,p,t) 

(3) 

0 = 

rt>(x(ti),p,U) 

(4) 

V = 

N(x,u,d, p) 

(5) 


7 Constrained Linearization 


Define the reference trajectory, x 0 
d = 0, p = 0, A£i C = 0. 

(f) by choosing: 

x 0 = F(x a , u 0 , 0, 0) 

(6) 

x o{tic) ~ X{c 

(7) 

0 = G(x o ,u o ,0,6,t) 

(8) 

o = il>(x e (tj),6,tf) 

(9) 

Vo - H(x o ,Go,0,d) 

(10) 


These equations of motion are integrated segment 
by segment to obtain a reference trajectory. The 
segment transition times are denoted t \ , fj, . . . , 
and the end of the trajectory by tj . 

The state variation is defined at a fixed point in 
time by the difference between the reference 
trajectory and the perturbed trajectory. 

bx(t) = x(t) — x 0 (t) (11) 

The equations of motion are linearized by 
evaluating the first variation, in the sense of the 
calculus of variations. 

6x — F x Sx 4 " F^Su 4 “ Fjd + Fpp ( 12 ) 

0 = G x 6x + Gu$u 4 " G&d 4 * Gpp ( 13 ) 

Since the constraint equations have been chosen to 
be well behaved, Equation (13) can be solved for 
the perturbed controls. 

6S = -GZ 1 ( G x 6x + Gj + Gppj (14) 

Substituting (14) into (12) and collecting terms: 


bx 


(F, - GZ l G t )6x 



+(F d - F u GZ l G d )d 



+(F P — F u GZ 1 G p )p 

bx 

= 

F x 6x + F d d+F p p (15) 

where : 


F x — F x Ft, G u Gg i etc. 


The output variation is found in a similar fashion. 
by = N x 6x + Ndd + Npi> (16) 

The above linearized system is defined at each 
point in time along the reference trajectory, with 
the equations changing discontinuously at the 
segment transitions. Additionally, the segment 
transition times will vary as the trajectory is 
perturbed. 

Define the state differential dx(f,) as the difference 
between x 0 (f«) and x(f. 4- eft,*), the state at the t’th 
reference transition at time U and the state at the 
perturbed transition at time 4* eft,. Since the 
equations change at time tj, the distinction will be 
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made below whether the values are before (t ~ ) or 
after (tf) the transition. 

dx(ti) = x(t,- + dti) - x 0 (ti) 

= x(t ~ ) + i(f ," )dti - x 0 (ti) 

= 6x(t-) + F(t-)dU (17) 

= 6x(it) + F(tt)dU (18) 

where, F(t) = F(x,p,d,t) = | G=0 . See Figure 2. 



The portion of the perturbed trajectory with a 
light line in Figure 2, between the reference 
transition time (t,*) and the actual transition time 
(<, + dti), is fictitious but is required later for the 
expressions for terminal error. This segment is 
required since the time differential is allowed to be 
arbitrary. This fictitious trajectory segment could 
be removed for the time histories of trajectory 
error, but since this is not the main purpose of the 
method, that is not described. 

Equations (15) and (16) represent a linear dynamic 
system with time varying coefficients that describe 
the response of the trajectory errors to the 
disturbances, with transitions between segments 
described by equations (17), (18) and (19). The 
following section will discuss various uses of this 
system to analyze the trajectory sensitivities. 

The method of Constrained Linearization can be 
formulated with any appropriate independent 
variable instead of time, such as path distance. 


8 Uses of the Linearized System 


8.1 Trajectory Error History 


Figure 2: Errors During Segment Transition 


The differential of the capture condition must be 
zero to ensure that it is satisfied at the perturbed 
transition time. This will require a differential in 
the segment transition time to satisfy the capture 
condition. 




— dx (t i+1 ) + — P +—dt i+l 
d^i dipi _ 

¥ M ‘ i+ i ) + F 



^6x(<r+i) + ^j±P+ipidti + 1 = 0 

, ; dipi =, dipi 

where, ipi = ~Q^F[t i+ ,) + — 


Solving for the differential in segment transition 
time: 



The linear time varying system that describes the 
trajectory errors can be integrated across each 
segment from the initial conditions at time t c to 
time tj . 


n(t) 

6x(t) = + 

1=1 

+ J* (f x 6x + Fid + F p p) dt (20) 

This can be integrated to determine the variation 
in the state history. The summation terms account 
for the change in the closed loop state equations 
across the segment transitions combined with the 
differential in the transition time. n(t) represents 
the number of transitions that have been 
encountered through time t , Substituting the 
definition of 6x from Equations (17) and (18): 

6x(t) = dx(t c )-F(t 0 )dt 0 
n(t) 

+DW> -*«■)]* 

i=i 
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+ J* (F z 6x + F d d + F p p) dt (21) 


Integrating by parts and collecting terms. 

**(*<+! ) = ) 


( 23 ) 


Note that the integration occurs from the reference 
initial time t 0 to the reference final time tf , with 
each segment being integrated from f, to since 
the state variations across the segment transitions 
due to the transition time differences are accounted 
for by the Fdti terms. (See Figure 2.) 

The differential in the final state can be found from 
Equation (17). 

dx(t s ) = 6z(tp + F(tf)dt f (22) 


8.2 Expressions for Terminal Error 


It is valuable to have expressions for the terminal 
error that do not require integrating the full linear 
time varying system for each set of disturbances 
and parameters chosen. This section describes an 
expression for the error at the end of a single 
segment based on the conditions at the beginning 
of the segment and the errors introduced during 
the segment. This method is then extended for 
multi-segment trajectories by accounting for the 
change in errors across the segment transitions. 


Evaluating Equation (20) from an initial time f, to 
a final time f j+i : 


«?(<r + i) = 6x(tf) 


+ Jt ( F t 6x + F d d + Fppj dt 


To obtain an expression for the final state variation 
i$ terms of only the disturbances and parameters, 
the state variation in the integrand of the above 
equation must be eliminated. This will be 
accomplished by adjoining the integrand with a 
Lagrange multiplier matrix A times the state error 
evolution equation. This does not change the 
integral, since the state equation is always zero, but 
is a trick that can be used to simplify the problem. 
The Lagrange multiplier will be chosen to eliminate 
6z in the integrand. 


- m +1 )6x(t- +l ) + A (tf)Sx(tt) 
+ f*’ +i [[(I + A)F X + A) 6x 

+ (/ + A)F d d + (/ + A)F p pj dt 

Choose: = 0 (24) 

A = -(/ + A)F, (25) 


A can be integrated backwards in time as a 
function of the reference trajectory alone - 
independent of the variations in the initial and final 
conditions, and the disturbances and parameters. 

A physical interpretation of A is the sensitivity of 
the final state to a variation in the current state 
derivative. A(<) = 


Substituting Equations (24) and (25) into (23), and 
gathering terms, the state variation at time f^. 1 is: 


**(<<+ 1 ) = °i 6£ ( t t )+*P 



where: 


a, = (I + Mtf)) 
bi(t) = (/ + A) Fj 



(26) 

(27) 

(28) 

(29) 


The variation in the final state at the final time is 
not the desired answer, since a differentia] in final 
state is required to satisfy the terminal condition. 
For a single segment trajectory, or the final 
segment of a multi-segment trajectory, the final 
state differential is found by extrapolating along 
the trajectory for the time differential (19). 


dx(t j+ i) 

dx(t i+ i) 


6x(t~ +1 ) + F(t~ +1 )dt i+1 


6x(t w ) = 
+ 
+ 


Sx(tf) 


J ‘ +l \F z 6x + F d d+ F p p 

e. 

A(F x Sx + F d d + Fpf— 65T ) J dt 

dx(t i+ i) 




*.• ^(*r+i) +e< p 


(30) 

(31) 

(32) 
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The matrix $j can be interpreted as taking the 
variation in state at the final time along the 
constrained equations of motion until the terminal 
condition is met. The matrix ©,• can be interpreted 
as the change in final state due to a perturbation of 
the termination condition tp. 


Substituting Equation (26) into (32), and 
simplifying the result. 


dx(t i+ i) = Ai 6x(tf ) + Qp 


+ / Bi(t)d{t) dt 

Jti 

(33) 

where: 


Ai = <t>i ai 

(34) 

Bi(t) = 

(35) 

Ci = <!>; c, + e, 

(36) 


For a segment transition, the errors are 
extrapolated from the end of the first segment to 
the perturbed transition time, then back to the 
reference transition time. Equating Expressions 
(17), (18) at time t,+i, substituting Equation (19), 
and solving for 6z(t++ i) yields: 


6*{tt+ i ) = *• + ©• P 

= / + (W + 
e. = (F(<t +1 )-F(<- +1 )) 

6x(tf +1 ) = Ai6x(tt) + Cip 


+ 



dt 


where, >1, , £,, C, are found from Equations (34) - 
(36). 


Extending the result to multiple segments: 
Ai = a, Ai-\ 

Bi(t) = (MO + «,*!- 0 

Ci = (Ci -I- CLiCi-l) + ©i 


conditions and parameters. B(f) can be examined 
as a function of time to observe when the trajectory 
is most sensitive to the disturbances, or integrated 
with a given disturbance history d(t) to obtain the 
total error due to a particular disturbance. 


8.3 Worst Case Disturbances 

The influence function B(t) indicates the influence 
of the disturbances at a point in time on the final 
states. This influence function can be used at each 
point in time to select a worst case disturbance to 
drive the final state error in a particular direction. 
For example, the disturbance vector may be 
magnitude limited and we are interested in what is 
the worst case error for any disturbance that meets 
this condition. The product of £(<) and the so 
chosen worst case disturbance d(f) can be 
integrated to obtain an estimate of the worst case 
errors for a given class of disturbances. 


8.4 Covariance Analysis 

The linear time varying system lends itself well to 
the analysis of random disturbances using 
covariance analysis. This analysis provides the 
standard deviation of the state error as a function 
of time given statistical properties of the 
disturbances. 

Bryson and Ho, [13] pages 342-345, describe a 
method of covariance analysis for a linear time 
varying system with random noise inputs. 

A method to account for the covariance 
propagation across segment transitions is under 
development. 


9 Example Application 


The expressions A , B(t) } and C can be interpreted 
as the influence functions of the initial conditions, 
disturbances, and parameters on the final state 
differential. The influence functions provide a 
quantification of the effects of the error sources on 
the final state differential. The influence functions 
A and C can be interpreted directly as the 
sensitivity of the final states to the initial 


An example is presented for a descent trajectory 
with a flight profile similar to a CTAS descent 
through an Approach Control Center airspace. The 
aircraft model used is the CTAS Boeing 757 model. 
The reference trajectory is define by the sequence 
of flight segments described in Table 1 below. Also 
shown in Table 1 are the segment transition 
conditions, the times the transitions occur, and the 
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changes to the transition times due to a drag 
perturbation. The reference trajectory is shown in 
Figure 3. 

The trajectory deviation due to a modeling error in 
drag of 5% is shown in Figure 4, as computed by 
the method described in Section 8.1. 

The sensitivity of the final states of the reference 
trajectory to the initial conditions and parametric 
error sources can be computed directly from the 
trajectory as described in Section 8.2. Table 2 
shows selected columns from the matrix A and 
Table 3 shows the matrix C. The columns of both 
tables have been scaled by reasonable magnitudes 
of the error sources for a better comparison of 
relative importance of trajectory error due to the 
error sources. Note that the reference trajectory 
ends in a level flight segment at specified airspeed 
and inertial position , so the Anal state errors are 
absorbed primarily in the final time. The idle 
thrust perturbation is between the nominal model 
of idle thrust and zero idle thrust at all flight 
conditions. 


10 Conclusions 


This paper develops a constrained linearization 
approach for the analysis of differences between 
computed and actual aircraft trajectories in air 
traffic control automation. Synthesized trajectories 
are needed in ATC automation to predict aircraft 
flight paths. Due to the modeling errors in 
trajectory computations, synthetic trajectories are 
different from actual trajectories. In this paper, 
modeling errors are described with uncertain 
parameters and uncertain functions. Pilot feedback 
control actions are expressed as equality 
constraints along a flight path. A constrained 
linearization approach is developed to analyze 
trajectory perturbations along a reference path 
caused by modeling errors. A numerical example is 
presented to illustrate the applications of the 
proposed techniques. The example demonstrates 
the efficiency and physical insight offered by the 
proposed approach. 
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Seg# 

Control Obj. 1 

Control Obj. 2 

Capture Cond. 

End Time 

Time Diff 

1 

h p = 35000 ft 

M = 0.78 

Xi = 70 nm 

560.8 sec 

0 sec 

2 

M = 0.78 

it = idle 

Vcas = 307 knt 

678.4 sec 

-5.6 sec 

3 

Vcas = 310 knt 

ir = idle 

hp = 10400 ft 

1094.0 sec 

-26.9 sec 

4 

h p = 10000 ft 

tt = idle 

Vcas = 250 knt 

1164.2 sec 

-30.8 sec 

5 

h p = 10000 ft 

Vcas = 250 knt 

Xi = 150 nm 

1312.3 sec 

-4-12.5 sec 


Table 1: Example Trajectory Description (A 757 Descent) 


Final 

State 

Error 

Initial State Error 

Vi 

1 knot 

Xi 

0.2 nmile 

Ml 

mg 

7000 lbs 

V, 

knots 

MM 

0 1 

0 

0 

7 a 

deg 

Bfl 

0 

0 

0 

a 

nmile 

■■ 

0 

0 

0 

mm 

feet 

mm 

0 

0 

0 


lbs 

B9 

2.80 

5.72 


Hi 

sec 


-1.60 

0 



Table 2: Terminal Sensitivity to Initial Conditions 


Final 

State 

Error 

Parametric Error Source 

Idle Thrust 
see text 

Drag 
+ 5% 

Wind 
5 knots 

VcASc 

2 knots 

M c 

0.005 

IS 

xtod 
2 nmiles 

Vt 

knot 

0 

0 

0 

2.27 

IKK? 

0.854 

0 

7 a 

deg 

0 

0 

0 

0 


0 

0 

X{ 

nmile 

0 

0 

0 

0 


0 

0 

hi 

feet 

0 

0 

0 

0 


200 

0 

mg 

lbs 

-27.6 


-33.8 

-10.3 


3.61 

15.9 

V 

sec 

2.0 


19.3 

-1.3 

ESI 

-1.04 

-8.9 


Table 3: Terminal Sensitivity to Parametric Errors 
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Figure 3: Boeing 757 Descent Reference Trajectory and Perturbed TVajectory 
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Figure 4: Trajectory Deviation Due to 5 % Drag Error 








